clear; clc;

load('data/video/escalator.mat');

data = double(data);
data = data/255;
[m, n, k] = size(data);
data = reshape(data, m*n, k);
% idx = ceil(rand(400, 1)*size(data, 2));
% data = data(:,idx);
% idx = ceil(rand(2000, 1)*size(data, 1));
% data = data(idx, :);

colMean = mean(data, 1);
data = bsxfun(@minus, data, colMean);

O = data + 0.15*randn(size(data));
colMean = mean(O, 1);
O = bsxfun(@minus, O, colMean);

PSNRin = psnr(O, data, 1);

traIdx = ceil(k*rand(ceil(k/2), 1));
tstIdx = ceil(k*rand(ceil(k/2), 1));

clear colMean m n k idx;

lambdamax = topksvd(O, 1, 5);
gridLambda = lambdamax*(3/4).^(10:16);
para.tol = 1e-10*numel(data);
para.tau = 1.01;
% % APG
% para.maxIter = 100;
% para.objstep = 1;
% para.regType = 1;
% for i = 1:length(gridLambda)
%     lambda = gridLambda(i);
%     mu = 0.5*lambda/sqrt(max(size(O)));
%     
%     [ X, Y ] = FastProxRPCA(O(:,traIdx), lambda, mu, 1e+8, 1e+8, para);
%     D = X + Y;
%     gridPSNR(i) = psnr(D, data(:,traIdx), 1);
%     clear X Y D;
%     
%     if(i > 1 && gridPSNR(i) < gridPSNR(i - 1))
%         break;
%     end
% end
% [~, lambda] = max(gridPSNR);
% lambda = gridLambda(lambda);
% mu = 0.5*lambda/sqrt(max(size(O)));
% 
% t = tic;
% [ X, Y, out(1) ] = APGRPCA(O(:,tstIdx), lambda, mu, para);
% Time(1) = toc(t);
% D = X + Y;
% PSNR(1) = psnr(D, data(:,tstIdx), 1);
% 
% clear X Y D;

for regType = 1:1
    
para.regType = regType;
para.maxIter = 50;
para.objstep = 1;

gridPSNR = zeros(size(gridLambda));
for i = 1:length(gridLambda)
    lambda = gridLambda(i);
    mu = lambda/sqrt(max(size(O)));
    
    switch(para.regType)
        case 1
            theta1 = lambda*1.25;
        case 2
            gridLambda(i) = 2*gridLambda(i);
            lambda = gridLambda(i);
            theta1 = sqrt(lambda);
        case 3
            theta1 = 10;
    end
    theta2 = 1e+8;
    
    [ X, Y ] = FastProxRPCA(O(:,traIdx), lambda, mu, theta1, theta2, para);
    D = X + Y;
    gridPSNR(i) = psnr(D, data(:,traIdx), 1);
    
    clear D X Y;
    
    if(i > 1 && gridPSNR(i) < 0.999*gridPSNR(i - 1))
        %break;
    end
end

[~, lambda] = max(gridPSNR);
lambda = gridLambda(lambda);

gridTheta = [4, 2, 1, 0.5, 0.25, 0.125];
gridPSNR = zeros(size(gridTheta));
for i = 1:length(gridTheta)
    theta1 = gridTheta(i)*lambda;
    theta2 = 1e+8;
    
    [ X, Y ] = FastProxRPCA(O(:,traIdx), lambda, mu, theta1, theta2, para);
    D = X + Y;
    gridPSNR(i) = psnr(D, data(:,traIdx), 1);
    
    if(i > 1 && gridPSNR(i) < gridPSNR(i - 1))
        break;
    end
end

[~, theta1] = max(gridPSNR);
theta1 = gridTheta(theta1);

t = tic;
[ X, Y, out(2,regType) ] = FastProxRPCA( O(:,tstIdx), lambda, mu, theta1, theta2, para);
Time(2,regType) = toc(t);
D = X + Y;
PSNR(2,regType) = psnr(D, data(:,tstIdx), 1);

clear X Y D;

%% GIST
t = tic;
[ X, Y, out(3,regType) ] = GISTRPCA( O(:,tstIdx), lambda, mu, theta1, theta2, para);
Time(3,regType) = toc(t);
D = X + Y;
PSNR(3,regType) = psnr(D, data(:,tstIdx), 1);

clear X Y D;

save('temp', 'out', 'PSNR', 'Time');

end

clear data O t traIdx tstIdx gridLambda gridMu i regType;

save('testVideo');
